---
title: "Replication File (Main Analyses): The Limits of Public Support for Fiscal Consolidation: Survey Evidence from Great Britain"
author: "Björn Bremer, Charlotte Cavaille, Lisanne de Blok, Catherine de Vries"
date: "February 3, 2026"
output:
  bookdown::html_document2:
    keep_tex: true
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE)
knitr::opts_chunk$set(message = FALSE)
knitr::opts_chunk$set(warning = FALSE)
```

```{r, include=FALSE}
rm(list = ls())
set.seed(1234567)

# List of packages needed
packages <- c("tidyverse",
              "haven",
              "bookdown",
              "stargazer",
              "gridExtra",
              "ggpubr",
              "margins",
              "coefplot",
              "modelsummary",
              "xtable",
              "estimatr",
              "kableExtra",
              "ltm", 
              "patchwork",
              "ggrepel"
)

# Function to check if packages are installed and install if needed
package_check <- function(pkg){
  if (!require(pkg, character.only = TRUE)) {
    install.packages(pkg, dependencies = TRUE)
    library(pkg, character.only = TRUE)
  }
}

# Apply function to each package
invisible(sapply(packages, package_check))

# Create directories to save data
if (!dir.exists("data")) {
  dir.create("data")
}

# Create directories to save figures
if (!dir.exists("figures")) {
  dir.create("figures")
}

```

# Load data from all waves

```{r, results = 'asis', message=F, echo=FALSE}
## Load all data

# load and clean debt data
debt_data <- read_csv("clean_general_gov_debt.csv") %>%
  dplyr::select(-1) %>%  # Drop first column (v1)
  rename(Debt = Value) %>%
  mutate(Debt = ifelse(Debt == "NA", NA, Debt),
         Debt = as.numeric(Debt)) %>%
  filter(Year > 1980)

# load and clean interest rates data
interest_data <- read_csv("clean_interest_rates.csv") %>%
  dplyr::select(-1) %>%  # Drop first column (v1)
  rename(Interest = value) %>%
  rename(Year = year) %>%
  mutate(Interest = ifelse(Interest == "NA", NA, Interest),
         Interest = as.numeric(Interest)) %>%
  filter(Year > 1980)

load("data/data_processed_w1.RData")
load("data/data_processed_w2.RData")
load("data/data_processed_w2_w3.RData")
```

# Preliminary analyes

```{r, results = 'asis', message=F, echo=FALSE}
# calculate Cronbach's alpha for items that are part of the debt aversion index 
w1_index <- subset(debt_w1, select = c(pol1_rec_w1, pol2_rec_w1, pol3_rec_w1))
cronbach.alpha(w1_index, na.rm = TRUE)
# result: 0.775
detach("package:ltm", unload=TRUE)
```

# Analyses for Main Paper

## Figure 1: Government debt and interest rates

```{r, results = 'asis', message=F, echo=FALSE}
# ============================================================================
# FIGURE 1A: Western/Southern European countries - Debt
# ============================================================================

# Create data with last point for labeling
last_points <- debt_data %>%
  filter(country %in% c("France", "United Kingdom", "Spain", "Portugal", "Austria", "Germany")) %>%
  group_by(country) %>%
  filter(Year == max(Year)) %>%
  ungroup()

fig1a <- debt_data %>%
  filter(country %in% c("France", "United Kingdom", "Spain", "Portugal", "Austria", "Germany")) %>%
  ggplot(aes(x = Year, y = Debt, color = country, linetype = country)) +
  geom_line(linewidth = 0.8) +
  geom_text_repel(data = last_points, 
                  aes(label = country), 
                  nudge_x = 1,
                  hjust = 0,
                  size = 2.5,
                  direction = "y",
                  segment.size = 0.2) +
  scale_color_manual(
    values = c("Austria" = "#000000",      # Black
               "France" = "#E69F00",       # Orange
               "Germany" = "#56B4E9",      # Sky blue
               "Portugal" = "#009E73",     # Bluish green
               "Spain" = "#CC79A7",        # Reddish purple
               "United Kingdom" = "#D55E00")  # Vermillion
  ) +
  scale_linetype_manual(
    values = c("Austria" = "solid", 
               "France" = "dashed", 
               "Germany" = "dotted", 
               "Portugal" = "longdash",
               "Spain" = "dotdash", 
               "United Kingdom" = "twodash")
  ) +
  scale_y_continuous(limits = c(0, 140), breaks = seq(0, 140, 20)) +
  scale_x_continuous(
    limits = c(1980, NA),
    breaks = seq(1990, 2020, 10),
    expand = expansion(mult = c(0.02, 0.15))
  ) +
  labs(x = "", y = "Public debt as a percent of GDP") +
  theme_minimal() +
  theme(
    legend.position = "none",
    axis.title.y = element_text(size = 9)
  )

# ============================================================================
# FIGURE 1B: Nordic countries - Debt
# ============================================================================

# Create data with last point for labeling
last_points_nordic <- debt_data %>%
  filter(country %in% c("Denmark", "Finland", "Norway", "Sweden")) %>%
  group_by(country) %>%
  filter(Year == max(Year)) %>%
  ungroup()

fig1b <- debt_data %>%
  filter(country %in% c("Denmark", "Finland", "Norway", "Sweden")) %>%
  ggplot(aes(x = Year, y = Debt, color = country, linetype = country)) +
  geom_line(linewidth = 0.8) +
  geom_text_repel(data = last_points_nordic, 
                  aes(label = country), 
                  nudge_x = 1,
                  hjust = 0,
                  size = 2.5,
                  direction = "y",
                  segment.size = 0.2) +
  scale_color_manual(
    values = c("Denmark" = "#E69F00",       # Orange
               "Finland" = "#56B4E9",       # Sky blue
               "Norway" = "#009E73",        # Bluish green
               "Sweden" = "#CC79A7")        # Reddish purple
  ) +
  scale_linetype_manual(
    values = c("Denmark" = "solid", 
               "Finland" = "dashed", 
               "Norway" = "dotted",
               "Sweden" = "longdash")
  ) +
  scale_y_continuous(limits = c(0, 140), breaks = seq(0, 140, 20)) +
  scale_x_continuous(
    limits = c(1980, NA),
    breaks = seq(1990, 2020, 10),
    expand = expansion(mult = c(0.02, 0.15))
  ) +
  labs(x = "", y = "") +
  theme_minimal() +
  theme(
    legend.position = "none",
    axis.title.y = element_text(size = 9)
  )

# ============================================================================
# FIGURE 1C: Western/Southern European countries - Interest Rates
# ============================================================================

# Filter for long-term interest rates only
lt_interest <- interest_data %>%
  filter(Measure == "Long-term interest rates")

# Create data with last point for labeling - Western/Southern European
last_points_interest <- lt_interest %>%
  filter(country %in% c("France", "United Kingdom", "Spain", "Portugal", "Austria", "Germany")) %>%
  group_by(country) %>%
  filter(Year == max(Year)) %>%
  ungroup()

fig1c <- lt_interest %>%
  filter(country %in% c("France", "United Kingdom", "Spain", "Portugal", "Austria", "Germany")) %>%
  ggplot(aes(x = Year, y = Interest, color = country, linetype = country)) +
  geom_line(linewidth = 0.8) +
  geom_text_repel(data = last_points_interest, 
                  aes(label = country), 
                  nudge_x = 1,
                  hjust = 0,
                  size = 2.5,
                  direction = "y",
                  segment.size = 0.2) +
  scale_color_manual(
    values = c("Austria" = "#000000",      # Black
               "France" = "#E69F00",       # Orange
               "Germany" = "#56B4E9",      # Sky blue
               "Portugal" = "#009E73",     # Bluish green
               "Spain" = "#CC79A7",        # Reddish purple
               "United Kingdom" = "#D55E00")  # Vermillion
  ) +
  scale_linetype_manual(
    values = c("Austria" = "solid", 
               "France" = "dashed", 
               "Germany" = "dotted", 
               "Portugal" = "longdash",
               "Spain" = "dotdash", 
               "United Kingdom" = "twodash")
  ) +
  scale_y_continuous(limits = c(0, 20), breaks = seq(0, 20, 2)) +
  scale_x_continuous(
    limits = c(1980, NA),
    breaks = seq(1990, 2020, 10),
    expand = expansion(mult = c(0.02, 0.15))
  ) +
  labs(x = "", y = "Long-term interest rates (10-year govt bonds)") +
  theme_minimal() +
  theme(
    legend.position = "none",
    axis.title.y = element_text(size = 9)
  )

# ============================================================================
# FIGURE 1D: Nordic countries - Interest Rates
# ============================================================================

# Create data with last point for labeling - Nordic countries
last_points_interest_nordic <- lt_interest %>%
  filter(country %in% c("Denmark", "Finland", "Norway", "Sweden")) %>%
  group_by(country) %>%
  filter(Year == max(Year)) %>%
  ungroup()

fig1d <- lt_interest %>%
  filter(country %in% c("Denmark", "Finland", "Norway", "Sweden")) %>%
  ggplot(aes(x = Year, y = Interest, color = country, linetype = country)) +
  geom_line(linewidth = 0.8) +
  geom_text_repel(data = last_points_interest_nordic, 
                  aes(label = country), 
                  nudge_x = 1,
                  hjust = 0,
                  size = 2.5,
                  direction = "y",
                  segment.size = 0.2) +
  scale_color_manual(
    values = c("Denmark" = "#E69F00",       # Orange
               "Finland" = "#56B4E9",       # Sky blue
               "Norway" = "#009E73",        # Bluish green
               "Sweden" = "#CC79A7")        # Reddish purple
  ) +
  scale_linetype_manual(
    values = c("Denmark" = "solid", 
               "Finland" = "dashed", 
               "Norway" = "dotted",
               "Sweden" = "longdash")
  ) +
  scale_y_continuous(limits = c(0, 20), breaks = seq(0, 20, 2)) +
  scale_x_continuous(
    limits = c(1980, NA),
    breaks = seq(1990, 2020, 10),
    expand = expansion(mult = c(0.02, 0.15))
  ) +
  labs(x = "", y = "") +
  theme_minimal() +
  theme(
    legend.position = "none",
    axis.title.y = element_text(size = 9)
  )

# ============================================================================
# COMBINE AND EXPORT
# ============================================================================

# Combine all four plots in 2x2 grid
combined_plot <- (fig1a | fig1b) / (fig1c | fig1d)

# Save as PDF
ggsave(
  filename = file.path("figures/Figure1.pdf"),
  plot = combined_plot,
  width = 18,
  height = 14,
units = c("cm"),
device = "pdf"
)

```

## Figure 2: Share of respondents who prioritize debt over competing issues

```{r, results = 'asis', message=F, echo=FALSE}
# subset to only two parties
debt_w1_rec <- debt_w1[which(debt_w1$party_rec_w1=="Labour" | debt_w1$party_rec_w1=="Conservatives"),]
debt_w2_w3_rec <- debt_w2_w3[which(debt_w2_w3$party_rec_w2=="Labour" | debt_w2_w3$party_rec_w2=="Conservatives"),]

# subset only to control group
debt_w1_rec <- debt_w1_rec[which(debt_w1_rec$treat_w1==1),]
debt_w2_w3_rec <- debt_w2_w3_rec[which(debt_w2_w3_rec$treat_w2!=2 & debt_w2_w3_rec$treat_w2!=3),]

# priority given to debt...
debt_CC_rec_w1_df <- datasummary_crosstab(debt_CC_rec_w1 ~ party_rec_w1,statistic = 1 ~ Percent("col") , data = debt_w1_rec, output = 'data.frame')
debt_immi_rec_w1_df <- datasummary_crosstab(debt_immi_rec_w1 ~ party_rec_w1,statistic = 1 ~ Percent("col") , data = debt_w1_rec, output = 'data.frame')
debt_nhs_rec_w1_df <- datasummary_crosstab(debt_nhs_rec_w1 ~ party_rec_w1,statistic = 1 ~ Percent("col") , data = debt_w1_rec, output = 'data.frame')

debt_CC_rec_w1_df <- debt_CC_rec_w1_df %>% rename(priority=debt_CC_rec_w1)
debt_CC_rec_w1_df$Topic <- "Climate change"
debt_immi_rec_w1_df <- debt_immi_rec_w1_df %>% rename(priority=debt_immi_rec_w1)
debt_immi_rec_w1_df$Topic <- "Immigration"
debt_nhs_rec_w1_df <- debt_nhs_rec_w1_df %>% rename(priority=debt_nhs_rec_w1)
debt_nhs_rec_w1_df$Topic <- "NHS"

df_priority <- rbind(debt_CC_rec_w1_df[,-2], debt_immi_rec_w1_df[,-2], debt_nhs_rec_w1_df[,-2])
df_priority <- df_priority %>% 
  pivot_longer(
    cols = `Labour`:`Conservatives`, 
    names_to = "party",
    values_to = "value")
df_priority <- df_priority %>% subset(priority!="All") %>% mutate_at(4, as.numeric)
df_priority$value <- df_priority$value/100
df_priority$priority <- ifelse(df_priority$priority==0, "Debt is not priority", ifelse(df_priority$priority==1, "Debt is priority", NA))

gg_prio_simple <- ggplot(df_priority[which(df_priority$priority == "Debt is priority"),], aes(x = party, y = value, fill = party)) +
  geom_bar(stat = "identity") +
  scale_fill_manual(values = c("grey40", "grey80")) +
  facet_grid(~ Topic, scales = "free") +
  labs(y = "Share", x = "Partisanship") +
  theme_minimal() +
  theme(legend.position = "none") +
  geom_text(aes(label = value), size = 2.5, color = "black", position = position_stack(vjust = 0.5)) +
  ylim(0, 0.61)
#gg_prio_simple

# how often given prio? 
debt_w1_rec$sum_debt_prio <- debt_w1_rec$debt_CC_rec_w1+ debt_w1_rec$debt_immi_rec_w1 + debt_w1_rec$debt_nhs_rec_w1
debt_w1_rec$sum_debt_prio_rec <- ifelse(debt_w1_rec$sum_debt_prio==3, "Debt over all 3", ifelse(debt_w1_rec$sum_debt_prio==2, "Debt over at least 2", ifelse(debt_w1_rec$sum_debt_prio==0, "Never prio debt", NA)))
sum_debt_prio_rec_w1_df <- datasummary_crosstab(sum_debt_prio_rec ~ party_rec_w1,statistic = 1 ~ Percent("col") , data = debt_w1_rec, output = 'data.frame')

sum_debt_prio_rec_w1_df <- sum_debt_prio_rec_w1_df %>% 
  pivot_longer(
    cols = `Labour`:`Conservatives`, 
    names_to = "party",
    values_to = "value")
sum_debt_prio_rec_w1_df <- sum_debt_prio_rec_w1_df %>% subset(sum_debt_prio_rec!="All") %>% mutate_at(4, as.numeric)
sum_debt_prio_rec_w1_df$value <- sum_debt_prio_rec_w1_df$value/100
sum_debt_prio_rec_w1_df<- sum_debt_prio_rec_w1_df[,-2]

gg_prio_comp <- ggplot(sum_debt_prio_rec_w1_df, aes(x = party, y = value, fill = party)) + 
      geom_bar(stat = "identity") +
      scale_fill_manual(values = c("grey40", "grey80")) +
      facet_grid(~ sum_debt_prio_rec, scales = "free") + 
      labs(y = "Share", x = "Partisanship") + 
      theme_minimal() + 
      theme(legend.position = "none") +
      geom_text(aes(label = value),size = 2.5, color="black", position = position_stack(vjust = 0.5)) + ylim(0,0.61)
#gg_prio_comp

prio_graph <- ggarrange(gg_prio_simple, gg_prio_comp, nrow = 1, ncol=2)
prio_graph

ggsave(prio_graph, width = 25, height = 15, units = c("cm"), file ="figures/Figure2.pdf")
```

## Figure 3: Average points or votes allocated to different fiscal policies in budgetary tasks

```{r, results = 'asis', message=F, echo=FALSE}
# plot results from the allocation exercise 

# survey 1
# full sample
allocation_w1_all <- datasummary(alloc_debt_w1 + alloc_cuts_edu_w1 + alloc_cuts_benefits_w1 + alloc_tax_income_w1 + alloc_tax_vat_w1 ~ mean*Arguments(na.rm=TRUE), data = debt_w1_rec, output = 'data.frame', fmt = 3)

allocation_w1_all$sample <- "Full"
colnames(allocation_w1_all)[1] <- "Variable"

# only labour
allocation_w1_L <- datasummary(alloc_debt_w1 + alloc_cuts_edu_w1 + alloc_cuts_benefits_w1 + alloc_tax_income_w1 + alloc_tax_vat_w1 ~ mean*Arguments(na.rm=TRUE), data = debt_w1_rec[which(debt_w1_rec$party_rec_w1 == "Labour"),], output = 'data.frame', fmt = 3)

allocation_w1_L$sample <- "Labour"
colnames(allocation_w1_L)[1] <- "Variable"

# only conservatives
allocation_w1_C <- datasummary(alloc_debt_w1 + alloc_cuts_edu_w1 + alloc_cuts_benefits_w1 + alloc_tax_income_w1 + alloc_tax_vat_w1 ~ mean*Arguments(na.rm=TRUE), data = debt_w1_rec[which(debt_w1_rec$party_rec_w1 == "Conservatives"),], output = 'data.frame', fmt = 3)

allocation_w1_C$sample <- "Conservatives"
colnames(allocation_w1_C)[1] <- "Variable"

# combine the results
allocation_w1 <- rbind(allocation_w1_all, allocation_w1_L, allocation_w1_C)
allocation_w1$mean <- (as.numeric(allocation_w1$mean))
allocation_w1$sample <- factor(allocation_w1$sample , levels = c("Conservatives", "Labour", "Full"))
allocation_w1$Variable <- factor(allocation_w1$Variable , levels = c("alloc_cuts_edu_w1", "alloc_cuts_benefits_w1","alloc_tax_vat_w1", "alloc_tax_income_w1","alloc_debt_w1"))

# plot a graph allocation exercise
graph_allocation_w1 <- ggplot(allocation_w1, aes(x = mean, y = Variable, fill = sample)) + 
      geom_bar(stat = "identity", position = "dodge") +
      scale_fill_manual(values = c("grey40", "grey60", "grey80"), breaks = c("Full", "Labour","Conservatives"), labels = c("Full", "Labour","Conservatives")) + 
  scale_y_discrete(labels = c("Cutting edu spending", "Cutting benefits", "Increasing VAT", "Increasing inc tax\nprogressively", "Increasing debt")) +
      theme_minimal() +
      labs(y = "", x = "Average points", fill = NULL) +
      ggtitle("NHS Funding Exercise") + 
      theme(plot.title = element_text(face = "bold")) +
      xlim(0, 37)
```

```{r, results = 'asis', message=F, echo=FALSE}
# plot results from the qvsr  

# full sample
qvsr_all <- datasummary(qvsr_cc + qvsr_nhs + qvsr_inflation + qvsr_immi + qvsr_debt ~ mean*Arguments(na.rm=TRUE), data = debt_w2_w3_qvsr, output = 'data.frame', fmt = 3)

qvsr_all$sample <- "Full"
colnames(qvsr_all)[1] <- "Variable"

# only labour
qvsr_L <- datasummary(qvsr_cc + qvsr_nhs + qvsr_inflation + qvsr_immi + qvsr_debt ~ mean*Arguments(na.rm=TRUE), data = debt_w2_w3_qvsr[which(debt_w2_w3_qvsr$party_rec_w2 == "Labour"),], output = 'data.frame', fmt = 3)

qvsr_L$sample <- "Labour"
colnames(qvsr_L)[1] <- "Variable"

# only conservatives
qvsr_C <- datasummary(qvsr_cc + qvsr_nhs + qvsr_inflation + qvsr_immi + qvsr_debt ~ mean*Arguments(na.rm=TRUE), data = debt_w2_w3_qvsr[which(debt_w2_w3_qvsr$party_rec_w2 == "Conservatives"),], output = 'data.frame', fmt = 3)

qvsr_C$sample <- "Conservatives"
colnames(qvsr_C)[1] <- "Variable"

# combine the results
qvsr_dist <- rbind(qvsr_all, qvsr_L, qvsr_C)
qvsr_dist$mean <- (as.numeric(qvsr_dist$mean))
qvsr_dist$sample <- factor(qvsr_dist$sample , levels = c("Conservatives", "Labour", "Full"))

qvsr_dist$Variable <- factor(qvsr_dist$Variable, levels = c("qvsr_inflation","qvsr_immi","qvsr_nhs","qvsr_cc","qvsr_debt"))

# plot a graph for qvsr
graph_qvsr <- ggplot(qvsr_dist, aes(x = mean, y = Variable, fill = sample)) + 
      geom_bar(stat = "identity", position = "dodge") +
      scale_fill_manual(values = c("grey40", "grey60", "grey80"), breaks = c("Full", "Labour","Conservatives"), labels = c("Full", "Labour","Conservatives")) + 
      scale_y_discrete(labels = c("Inflation support","Immigration","NHS", "Climate change", "Debt")) +
      theme_minimal() +
      labs(y = "", x = "Average votes", fill = NULL) +
      ggtitle("QVSR") + 
      theme(plot.title = element_text(face = "bold")) +
      xlim(0, 3)
```

```{r, results = 'asis', message=F, echo=FALSE}
# combine the graphs
graph_exercises <- ggarrange(graph_allocation_w1, graph_qvsr, 
                         common.legend = T, legend = "bottom", nrow = 1, ncol=2)

# save the combined graph
ggsave(graph_exercises, width = 18, height = 12, units = c("cm"), file ="figures/Figure3.pdf")
```

## Figure 5: Average treatment effects concern, debt aversion index, and policy priorities

```{r, include = FALSE}
# analyze effect of treatment on debt concern
dd_debt_3 <- lm(debt_post_z_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_debt_2 <- lm(debt_post_z_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_debt_1 <- lm(debt_post_z_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine observations from all surveys in one graph
ate_debt_all <- multiplot(dd_debt_3, dd_debt_2, dd_debt_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Concern", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.87, 0.87), breaks = seq(-0.8, 0.8, 0.2)) +
      scale_color_grey(name = "Survey", breaks = c("dd_debt_3", "dd_debt_2", "dd_debt_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_debt_3", "dd_debt_2", "dd_debt_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_debt_3", "dd_debt_2", "dd_debt_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2")) + theme(plot.title = element_text(face = "bold"))
```

```{r, include = FALSE}
# analyze effect of treatment on debt index
dd_index_3 <- lm(debt_index_z_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_index_2 <- lm(debt_index_z_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_index_1 <- lm(debt_index_z_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine results from all waves in one graph
ate_index_all <- multiplot(dd_index_3, dd_index_2, dd_index_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Debt Aversion Index", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.87, 0.87), breaks = seq(-0.8, 0.8, 0.2)) +
      scale_color_grey(name = "Survey", breaks = c("dd_index_3", "dd_index_2", "dd_index_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_index_3", "dd_index_2", "dd_index_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_index_3", "dd_index_2", "dd_index_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"))  + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include = FALSE}
# analyze effect of treatment on debt vs. NHS
dd_nhs_3 <- lm(debt_nhs_rec_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_nhs_2 <- lm(debt_nhs_rec_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_nhs_1 <- lm(debt_nhs_rec_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine results from all waves in one graph
ate_nhs_all <- multiplot(dd_nhs_3, dd_nhs_2, dd_nhs_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Debt vs. NHS", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.2, 0.2), breaks = seq(-0.2, 0.2, 0.1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_nhs_3", "dd_nhs_2", "dd_nhs_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_nhs_3", "dd_nhs_2", "dd_nhs_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_nhs_3", "dd_nhs_2", "dd_nhs_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"))  + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include = FALSE}
# analyze effect of treatment on debt vs. climate change
dd_cc_3 <- lm(debt_CC_rec_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_cc_2 <- lm(debt_CC_rec_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_cc_1 <- lm(debt_CC_rec_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine results from all waves in one graph
ate_cc_all <- multiplot(dd_cc_3, dd_cc_2, dd_cc_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Debt vs. Climate Change", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.2, 0.2), breaks = seq(-0.2, 0.2, 0.1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_cc_3", "dd_cc_2", "dd_cc_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_cc_3", "dd_cc_2", "dd_cc_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_cc_3", "dd_cc_2", "dd_cc_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"))  + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include = FALSE}
# analyze effect of treatment on debt vs. immigration
dd_immi_3 <- lm(debt_immi_rec_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_immi_2 <- lm(debt_immi_rec_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_immi_1 <- lm(debt_immi_rec_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine results from all waves in one graph
ate_immi_all <- multiplot(dd_immi_3, dd_immi_2, dd_immi_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Debt vs. Immigration", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.2, 0.2), breaks = seq(-0.2, 0.2, 0.1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_immi_3", "dd_immi_2", "dd_immi_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_immi_3", "dd_immi_2", "dd_immi_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.2, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_immi_3", "dd_immi_2", "dd_immi_1"), labels = c("Survey 1", "Survey 2, wave 1", "Survey 2, wave 2"))  + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include = FALSE}
# analyze effect of treatment on debt vs. inflation
dd_infl_to_2 <- lm(debt_infl_rec_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

dd_infl_to_1 <- lm(debt_infl_rec_w3 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3)

# combine results from all waves in one graph
ate_infl_all <- multiplot(dd_infl_to_2, dd_infl_to_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Debt vs. Inflation Support", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.2, 0.2), breaks = seq(-0.2, 0.2, 0.1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_infl_to_2", "dd_infl_to_1"), labels = c("Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.4, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_infl_to_2", "dd_infl_to_1"), labels = c("Survey 2, wave 1", "Survey 2, wave 2"), start = 0.6, end = 0.4, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_infl_to_2", "dd_infl_to_1"), labels = c("Survey 2, wave 1", "Survey 2, wave 2"))  + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, results = 'asis', message=F, echo=FALSE}
# combine all graphs
ate_graph <- ggarrange(ate_debt_all, ate_index_all + 
                         theme(axis.text.y = element_blank(),
                               axis.title.y = element_blank()),
                       ate_cc_all, ate_immi_all + 
                             theme(axis.text.y = element_blank(),
                                   axis.title.y = element_blank()), 
                      ate_nhs_all, ate_infl_all + 
                             theme(axis.text.y = element_blank(),
                                   axis.title.y = element_blank()), 
                             common.legend = T, legend = "bottom", nrow = 3, ncol=2,
                             widths = c(1, 0.75))
ate_graph

# save the combined graph
ggsave(ate_graph, width = 18, height = 20, units = c("cm"), file ="figures/Figure5.pdf")
```

## Figure 6a: Average treatment effects on budgetary tasks

```{r, include = FALSE}
# analyse effect of treatment on reduction of debt in the NHS funding exercise 

# run regresion models
dd_alloc_3 <- lm(alloc_debt_z_w1 ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)

dd_alloc_2 <- lm(alloc_debt_z_w2 ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine results from all waves in one graph
ate_alloc_all <- multiplot(dd_alloc_3, dd_alloc_2, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "NHS Funding Exercise", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.515, 0.515), breaks = seq(-0.5, 0.5, 0.1), labels = round(seq(-0.5, 0.5, 0.1), 1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_alloc_3", "dd_alloc_2"), labels = c("Survey 1", "Survey 2"), start = 0.6, end = 0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_alloc_3", "dd_alloc_2"), labels = c("Survey 1", "Survey 2"), start = 0.6, end = 0, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_alloc_3", "dd_alloc_2"), labels = c("Survey 1", "Survey 2"))  + 
      theme(plot.title = element_text(face = "bold"))
ate_alloc_all
```

```{r, include = FALSE}
# analyse effect of treatment on debt in qvsr

# run regresion models
dd_qvsr_debt <- lm(qvsr_debt_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2_w3_qvsr)

# plot the results in a graph
ate_alloc_qvsr <- coefplot(dd_qvsr_debt,intercept = FALSE, coefficients = c("treat2", "treat3"), title = "QVSR", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.515, 0.515), breaks = seq(-0.5, 0.5, 0.1), labels = round(seq(-0.5, 0.5, 0.1), 1)) +
      scale_color_grey(name = "Survey", breaks = c("dd_qvsr_debt"), labels = c("Survey 2"), start = 0.6, end = 0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey",  breaks = c("dd_qvsr_debt"), labels = c("Survey 2"), start = 0.6, end = 0, aesthetics = "colour") +
      scale_shape_discrete(name = "Survey",  breaks = c("dd_qvsr_debt"), labels = c("Survey 2"))  + 
      theme(plot.title = element_text(face = "bold"))

ate_alloc_qvsr
```

```{r, results = 'asis', message=F, echo=FALSE}
# combine all graphs
ate_graph_alloc <- ggarrange(ate_alloc_all, ate_alloc_qvsr + 
                         theme(axis.text.y = element_blank(),
                               axis.title.y = element_blank()),
                             common.legend = T, legend = "bottom", nrow = 1, ncol=2,
                             widths = c(1, 0.75))
ate_graph_alloc

# save the combined graph
ggsave(ate_graph_alloc, width = 18, height = 8, units = c("cm"), file ="figures/Figure6a.pdf")
```

## Figure 6b: Average treatment effects on letter writing behavior

```{r}
# analyze effect of treatment on debt index
dd3_wrote <- lm(wroteF ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)
dd2_prioDebt <- lm(prioDebt ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)
dd1_prioSupp <- lm(prioSupp ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine results in one graph (order matters: wrote, prioDebt, prioSupp)
ate_letter_all <- multiplot(dd1_prioSupp, dd2_prioDebt, dd3_wrote,
                            intercept = FALSE, 
                            coefficients = c("treat2", "treat3"), 
                            title = "Letter Writing Outcomes", 
                            ylab = "Treatment", 
                            xlab = "Average treatment effect", 
                            pointSize = 2, 
                            innerCI = 0, 
                            newNames = c("treat2" = "Sustainable", 
                                        "treat3" = "Unsustainable"), 
                            dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(limits = c(-0.2, 0.2), breaks = seq(-0.2, 0.2, 0.1), labels = round(seq(-0.2, 0.2, 0.1), 1)) +
  scale_color_grey(name = "Outcome", 
                   breaks = c("dd3_wrote", "dd2_prioDebt", "dd1_prioSupp"), 
                   labels = c("Wrote Letter", "Prio: Debt", "Prio: Support"), 
                   start = 0.6, end = 0.2, aesthetics = "colour") +  
  scale_fill_grey(name = "Outcome",  
                   breaks = c("dd3_wrote", "dd2_prioDebt", "dd1_prioSupp"), 
                  labels = c("Wrote Letter", "Prio: Debt", "Prio: Support"), 
                  start = 0.6, end = 0.2, aesthetics = "colour") +
  scale_shape_discrete(name = "Outcome",  
                   breaks = c("dd3_wrote", "dd2_prioDebt", "dd1_prioSupp"), 
                       labels = c("Wrote Letter", "Prio: Debt", "Prio: Support")) + 
  theme(plot.title = element_text(face = "bold"),
        legend.position = "right") +
  guides(color = guide_legend(nrow = 3),
         shape = guide_legend(nrow = 3),
         fill = guide_legend(nrow = 3))

# Save plot
ggsave("figures/Figure6b.pdf", 
       ate_letter_all, 
       width = 14, 
       height = 7,
      units = c("cm"))
```


## Figure 7: Average treatment effects on expectations
```{r, include=FALSE}
## average treatment effect on expectations for nhs (quality)
debt_w1$exp_nhs_qual_w1 <- as.numeric(debt_w1$exp_nhs_qual_w1)
debt_w2$exp_nhs_qual_w2 <- as.numeric(debt_w2$exp_nhs_qual_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_nhs_qual_w1_z = as.vector(scale(exp_nhs_qual_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_nhs_qual_w2_z = as.vector(scale(exp_nhs_qual_w2)))

dd_nhs_qual_2 <- lm(exp_nhs_qual_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_nhs_qual_1 <- lm(exp_nhs_qual_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

ate_nhs_qual_all <- multiplot(dd_nhs_qual_2, dd_nhs_qual_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "NHS (Quality)", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3, 0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_nhs_qual_2", "dd_nhs_qual_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_nhs_qual_2", "dd_nhs_qual_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_nhs_qual_2", "dd_nhs_qual_1"), labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include=FALSE}
## average treatment effects on expectations for the nhs (costs)
debt_w1$exp_nhs_costs_w1 <- as.numeric(debt_w1$exp_nhs_costs_w1)
debt_w2$exp_nhs_costs_w2 <- as.numeric(debt_w2$exp_nhs_costs_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_nhs_costs_w1_z = as.vector(scale(exp_nhs_costs_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_nhs_costs_w2_z = as.vector(scale(exp_nhs_costs_w2)))

dd_nhs_costs_2 <- lm(exp_nhs_costs_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_nhs_costs_1 <- lm(exp_nhs_costs_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

ate_nhs_costs_all <- multiplot(dd_nhs_costs_2, dd_nhs_costs_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "NHS (Costs)", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3, 0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_nhs_costs_2", "dd_nhs_costs_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_nhs_costs_2", "dd_nhs_costs_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_nhs_costs_2", "dd_nhs_costs_1"),labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include=FALSE}
## average treatment effect on growth expectations
debt_w1$exp_gdp_w1 <- as.numeric(debt_w1$exp_gdp_w1)
debt_w2$exp_gdp_w2 <- as.numeric(debt_w2$exp_gdp_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_gdp_w1_z = as.vector(scale(exp_gdp_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_gdp_w2_z = as.vector(scale(exp_gdp_w2)))

dd_gdp_2 <- lm(exp_gdp_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_gdp_1 <- lm(exp_gdp_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine at in one graph
ate_gdp_all <- multiplot(dd_gdp_2, dd_gdp_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "GDP", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3, 0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_gdp_2", "dd_gdp_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.2, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_gdp_2", "dd_gdp_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.2, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_gdp_2", "dd_gdp_1"), labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include=FALSE}
## average treatment effect on inflation expectations
debt_w1$exp_infl_w1 <- as.numeric(debt_w1$exp_infl_w1)
debt_w2$exp_infl_w2 <- as.numeric(debt_w2$exp_infl_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_infl_w1_z = as.vector(scale(exp_infl_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_infl_w2_z = as.vector(scale(exp_infl_w2)))

dd_infl_2 <- lm(exp_infl_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_infl_1 <- lm(exp_infl_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine ates in one graph
ate_infl_all <- multiplot(dd_infl_2, dd_infl_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Inflation", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3, 0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_infl_2", "dd_infl_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_infl_2", "dd_infl_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_infl_2", "dd_infl_1"), labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include=FALSE}
## average treatment effect on tax expectation 
debt_w1$exp_tax_w1 <- as.numeric(debt_w1$exp_tax_w1)
debt_w2$exp_tax_w2 <- as.numeric(debt_w2$exp_tax_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_tax_w1_z = as.vector(scale(exp_tax_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_tax_w2_z = as.vector(scale(exp_tax_w2)))

dd_tax_2 <- lm(exp_tax_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_tax_1 <- lm(exp_tax_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine ates in one graph
ate_tax_all <- multiplot(dd_tax_2, dd_tax_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Taxes", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3,0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_tax_2", "dd_tax_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_tax_2", "dd_tax_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.0, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_tax_2", "dd_tax_1"),labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, include=FALSE}
## average treatment effect on pension expectations
debt_w1$exp_pen_w1 <- as.numeric(debt_w1$exp_pen_w1)
debt_w2$exp_pen_w2 <- as.numeric(debt_w2$exp_pen_w2)

# standardize variables
debt_w1 <- debt_w1 %>% 
  mutate(exp_pen_w1_z = as.vector(scale(exp_pen_w1)))
debt_w2 <- debt_w2 %>% 
  mutate(exp_pen_w2_z = as.vector(scale(exp_pen_w2)))

dd_pen_2 <- lm(exp_pen_w1_z ~ treat + as.factor(debt_pre_w1) + as.factor(block_w1), data = debt_w1)
dd_pen_1 <- lm(exp_pen_w2_z ~ treat + as.factor(debt_pre_w2) + as.factor(block_w2), data = debt_w2)

# combine ates in one graph
ate_pen_all <- multiplot(dd_pen_2, dd_pen_1, intercept = FALSE, coefficients = c("treat2", "treat3"), title = "Pensions", ylab = "Treatment", xlab = "Average treatment effect", pointSize = 2, innerCI=0, newNames = c("treat2" = "Sustainable", "treat3" = "Unsustainable"), dodgeHeight = 0.85) + 
      theme_bw() + scale_x_continuous(breaks = seq(from = -0.25, to = 0.25, by = 0.1), limits = c(-0.3, 0.3)) +
      scale_color_grey(name = "Survey", breaks = c("dd_pen_2", "dd_pen_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.1, aesthetics = "colour") +  
      scale_fill_grey(name = "Survey", breaks = c("dd_pen_2", "dd_pen_1"), labels = c("Survey 1", "Survey 2"), start = 0.4, end = 0.1, aesthetics = "colour") +
        scale_shape_discrete(name = "Survey", breaks = c("dd_pen_2", "dd_pen_1"), labels = c("Survey 1", "Survey 2")) + 
      theme(plot.title = element_text(face = "bold"))
```

```{r, results = 'asis', message=F, echo=FALSE}
# save combined plot
exp_all <- ggarrange(ate_tax_all, ate_pen_all +
                           theme(axis.text.y = element_blank(),
                                 axis.title.y = element_blank()),
                     ate_nhs_qual_all, 
                     ate_nhs_costs_all + 
                           theme(axis.text.y = element_blank(),
                                 axis.title.y = element_blank()),
                     ate_gdp_all, 
                     ate_infl_all +
                           theme(axis.text.y = element_blank(),
                                 axis.title.y = element_blank()),
                     common.legend = T, legend = "bottom", nrow = 3, ncol=2,
                     widths = c(1, 0.75))
exp_all

ggsave(exp_all, width = 18, height = 20, units = c("cm"), file ="figures/Figure7.pdf")
```

## Figure 8: HTEs by partisanship

```{r, results = 'asis', message=F, echo=FALSE}
# analyze effect of treatment on concern

# analyze effect of treatment on debt aversion index in survey 2, wave 2 (follow-up)
dd_index_part_w3 <- lm(debt_index_z_w3 ~ treat*party_rec_w2 + as.factor(debt_pre_w2), data = debt_w2_w3)
dd_index_part_ma_w3 <- margins(dd_index_part_w3, variables = "treat", at = list(party_rec_w2 = c("Labour", "Conservatives")))
ame_dd_index_part_w3 <- as.data.frame(summary(dd_index_part_ma_w3))
ame_dd_index_part_w3$party_rec_w2 <- as.character(ame_dd_index_part_w3$party_rec_w2)
ame_dd_index_part_w3$AME <- as.numeric(ame_dd_index_part_w3$AME)
ame_dd_index_part_w3$SE <- as.numeric(ame_dd_index_part_w3$SE)
ame_dd_index_part_w3$Treatment <- ame_dd_index_part_w3$factor
ame_dd_index_part_w3$upper95 <- ame_dd_index_part_w3$AME + (ame_dd_index_part_w3$SE*1.96)
ame_dd_index_part_w3$lower95 <- ame_dd_index_part_w3$AME - (ame_dd_index_part_w3$SE*1.96)

ame_dd_index_part_w3$Treatment[ame_dd_index_part_w3$Treatment=="treat2"] <- "Sustainable"
ame_dd_index_part_w3$Treatment[ame_dd_index_part_w3$Treatment=="treat3"] <- "Unsustainable"

ame_dd_index_part_w3$party_rec_w2 <- factor(ame_dd_index_part_w3$party_rec_w2, levels = c("Labour", "Conservatives"))
 
gg_ame_dd_index_part_w3 <- ggplot(ame_dd_index_part_w3, aes(colour = party_rec_w2)) + 
  geom_errorbar(aes(y = Treatment, xmin=lower95, xmax=upper95), width=0, position=position_dodge(width=0.5)) +
  geom_point(aes(y = Treatment, x=AME, shape = party_rec_w2), position = position_dodge(width=0.5), size=2) +
  theme_bw() + 
   scale_color_grey(name = "", breaks = c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"), start = 0.6, end = 0, aesthetics = "colour") +  
      scale_fill_grey(name = "",  breaks = c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"), start = 0.6, end = 0, aesthetics = "colour") +
      scale_shape_discrete(name = "",  breaks =c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"))  + 
  geom_vline(xintercept=0, color="grey", linetype = "dashed") +
      scale_x_continuous(limits = c(-0.65, 0.65), breaks = seq(-0.5, 0.5, 0.2)) +
      xlab("Average treatment effect") +
      theme(legend.position="bottom") +
      theme(legend.position="bottom") + ggtitle("Debt Aversion Index") + 
      theme(plot.title = element_text(face = "bold"))

# analyze effect of treatment on debt in qvsr in survey 2, wave 2 (follow-up)
dd_qvsr_debt_part <- lm(qvsr_debt_z ~ treat*party_rec_w2 + as.factor(debt_pre_w2), data = debt_w2_w3_qvsr)
dd_qvsr_debt_part_ma <- margins(dd_qvsr_debt_part, variables = "treat", at = list(party_rec_w2 = c("Labour", "Conservatives")))
ame_qvsr_debt_part <- as.data.frame(summary(dd_qvsr_debt_part_ma))
ame_qvsr_debt_part$party_rec_w2 <- as.character(ame_qvsr_debt_part$party_rec_w2)
#ame_dd_alloc_part_w2 <- rbind(ame_dd_alloc_part_w2, ate_row_tr2, ate_row_tr3)
ame_qvsr_debt_part$AME <- as.numeric(ame_qvsr_debt_part$AME)
ame_qvsr_debt_part$SE <- as.numeric(ame_qvsr_debt_part$SE)
ame_qvsr_debt_part$Treatment <- ame_qvsr_debt_part$factor
ame_qvsr_debt_part$upper95 <- ame_qvsr_debt_part$AME + (ame_qvsr_debt_part$SE*1.96)
ame_qvsr_debt_part$lower95 <- ame_qvsr_debt_part$AME - (ame_qvsr_debt_part$SE*1.96)

ame_qvsr_debt_part$Treatment[ame_qvsr_debt_part$Treatment=="treat2"] <- "Sustainable"
ame_qvsr_debt_part$Treatment[ame_qvsr_debt_part$Treatment=="treat3"] <- "Unsustainable"

ame_qvsr_debt_part$party_rec_w2 <- factor(ame_qvsr_debt_part$party_rec_w2, levels = c("Labour", "Conservatives"))

gg_ame_dd_qvsr_debt_part <- ggplot(ame_qvsr_debt_part, aes(colour = party_rec_w2)) + 
  geom_errorbar(aes(y = Treatment, xmin=lower95, xmax=upper95), width=0, position=position_dodge(width=0.5)) + 
  geom_point(aes(y = Treatment, x=AME, shape = party_rec_w2), position = position_dodge(width=0.5), size=2) +
  theme_bw() +
  scale_color_grey(name = "", breaks = c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"), start = 0.6, end = 0, aesthetics = "colour") +  
      scale_fill_grey(name = "",  breaks = c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"), start = 0.6, end = 0, aesthetics = "colour") +
      scale_shape_discrete(name = "",  breaks =c("Conservatives", "Labour"), labels = c("Conservatives", "Labour"))  + 
  geom_vline(xintercept=0, color="grey", linetype = "dashed") +
        scale_x_continuous(limits = c(-0.65, 0.65), breaks = seq(-0.5, 0.5, 0.2)) +
      xlab("Average treatment effect") +
      theme(legend.position="bottom") +
      theme(legend.position="bottom") + ggtitle("QVSR") + 
      theme(plot.title = element_text(face = "bold"))

# combine the figures
hte_partisan <- ggarrange(gg_ame_dd_index_part_w3, 
                              gg_ame_dd_qvsr_debt_part + 
                             theme(axis.text.y = element_blank(),
                                   axis.title.y = element_blank()),
                             nrow=1, ncol=2, common.legend = T,  legend="bottom",
                             widths = c(1, 0.75))
hte_partisan

ggsave(hte_partisan, width = 18, height = 8, units = c("cm"), file ="figures/Figure8.pdf")
```

# Letter Writing Exercise (results mentioned in main text)

```{r, results = 'asis', message=F, echo=FALSE}
# create variable whether people wrote something sensible 
debt_w2 <- debt_w2 %>%
  mutate(
    junk_not_useful = as.numeric(junk_not_useful),
    wrote_S = as.numeric(wrote_S),
    # Replace NA with 0 in junk_not_useful
    junk_not_useful = replace_na(junk_not_useful, 0)
  )

debt_w2 <- debt_w2 %>%
  mutate(wrote_S2 = ifelse(junk_not_useful == 1 & wrote_S == 1, 0, wrote_S)) 

# Cross-tabulation
table(debt_w2$wrote_S2, debt_w2$wroteF) # check that wrote_S2 is the same as wroteF - continue using wroteF

# Numbers mentioned in main text
mean(debt_w2$prioDebt, na.rm = TRUE)
mean(debt_w2$prioSupp, na.rm = TRUE)
mean(debt_w2$wrote_S2, na.rm = TRUE)
mean(debt_w2$junk_not_useful, na.rm = TRUE)

# Conditional on wroteF == 1: prioDebt
debt_w2_wrote <- debt_w2 %>% filter(wroteF == 1) 

mean(debt_w2_wrote$prioDebt, na.rm = TRUE)
debt_w2_wrote %>% filter(party_rec_w2 == "Labour") %>% summarise(mean_priod = mean(prioDebt, na.rm = TRUE))
debt_w2_wrote %>% filter(party_rec_w2 == "Conservatives") %>% summarise(mean_priod = mean(prioDebt, na.rm = TRUE))

# Conditional on wroteF == 1: prioSupp
mean(debt_w2_wrote$prioSupp, na.rm = TRUE)
debt_w2_wrote %>% filter(party_rec_w2 == "Labour") %>% summarise(mean_help = mean(prioSupp, na.rm = TRUE))
debt_w2_wrote %>% filter(party_rec_w2 == "Conservatives") %>% summarise(mean_help = mean(prioSupp, na.rm = TRUE))
```
